#ifndef _MOE_SBXCROSSOVER
#define _MOE_SBXCROSSOVER

#include <iostream>

using namespace std;

#include "engine/Random.hh"
#include "problems/RealArraySolution.hh"
#include "engine/Crossover.hh"

/**This crossover operator it's for apply to a RealArraySolution
  *It's an implementation of the Deb's SBX Crossover, see Simulated Binary Crossover For Continous Search Space,
  *by Deb and Ram Bhusan Agrawal November 1994. That paper describes a crossover operator that emulates the behavior of
  *a single point crossover over a real codification.
  */
class SBXCrossover: public Crossover {
private:
    double index; //Index of crossover
    static const double EPS = 1.0e-14;

public:

    SBXCrossover(double index = 20.0) : Crossover() {
        this->index = index;
    }

    double getIndex() {
        return index;
    }

    void setIndex(double index) {
        this->index = index;
    }

    virtual ~SBXCrossover() {}

    int matesExpected() {
        return 2;
    }

    ArrayList<Solution> *cross(ArrayList<Solution> *mates) {

        RealArraySolution *mate1 = (RealArraySolution *)mates->get(0);
        RealArraySolution *mate2 = (RealArraySolution *)mates->get(1);

        RealArraySolution *child0 = new RealArraySolution();
        RealArraySolution *child1 = new RealArraySolution();

        for (int i = 0; i < RealArraySolution::getSize(); i++) {
            double value1;
            double value2;

            combines(mate1->get(i), mate2->get(i), value1, value2, i);

            child0->set(i, value1);
            child1->set(i, value2);
        }

        ArrayList<Solution> *childs = new ArrayList<Solution>(false);
        childs->add(child0);
        childs->add(child1);

        if (RealArraySolution::getLowerBounds() == NULL || RealArraySolution::getHigherBounds() == NULL) {
            cout << endl << "RealArraySolution higher or lower bounds not defined, please set and try again, aborting!"
                 << endl;
            exit(1);
        }

        return childs;
    }

private:

    void combines(double value1, double value2, double &valueChild1, double &valueChild2, int idx) {

        if (Random::getInstance()->flip(0.5) && fabs(value1 - value2) > EPS) {
            double y1 = value1;
            double y2 = value2;

            if (y1 > y2) {
                double temp = y1;
                y1 = y2;
                y2 = temp;
            }

            //note that lower and higher bound must be not NULL, it's checked on SBXCrossover::cross
            double yl = RealArraySolution::getLowerBounds()[idx];
            double yh = RealArraySolution::getHigherBounds()[idx];
            double random = Random::getInstance()->getReal(0.0, 1.0);
            double beta = 1.0 + (2.0 * (y1 - yl) / (y2 - y1));
            double alpha = 2.0 - pow(beta, -(index + 1.0));

            double betaq = 0.0;

            if (random <= (1.0 / alpha))
                betaq = pow((random * alpha), (1.0 / (index + 1.0)));
            else
                betaq = pow((1.0 / (2.0 - random * alpha)), (1.0 / (index + 1.0)));

            valueChild1 = 0.5 * ((y1 + y2) - betaq * (y2 -y1));

            beta = 1.0 + (2.0 * (yh - y2) / (y2 - y1));
            alpha = 2.0 - pow(beta, -(index + 1.0));

            if (random <= (1.0 / alpha))
                betaq = pow((random * alpha), (1.0 / (index + 1.0)));
            else
                betaq = pow((1.0 / (2.0 - random * alpha)), (1.0 / (index + 1.0)));

            valueChild2 = 0.5 * ((y1 + y2) + betaq * (y2 - y1));

            if (Random::getInstance()->flip(0.5)) {
                double temp = valueChild1;
                valueChild1 = valueChild2;
                valueChild2 = temp;
            }

        } else {

            valueChild2 = value2;
            valueChild1 = value1;
        }
    }
};


#endif
